Crystallization in a sheared colloidal suspension 
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We study numerically the crystallization process in a supersaturated suspension of repulsive col- 
loidal particles driven by simple shear flow. The effect of the shear flow on crystallization is two-fold: 
while it suppresses the initial nucleation, once a large enough critical nucleus has formed its growth 
is enhanced by the shear flow. Combining both effects implies an optimal strain rate at which the 
overall crystallization rate has a maximum. To gain insight into the underlying mechanisms, we 
employ a discrete state model describing the transitions between the local structural configurations 
around single particles. We observe a time-scale separation between these transitions and the overall 
progress of the crystallization allowing for an effective Markovian description. By using this model, 
we demonstrate that the suppression of nucleation is due to the inhibition of a pre-structured liquid. 
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I. INTRODUCTION 

The freezing of a disordered colloidal or nano suspen- 
sion into a crystalline state with long-range order is a 
first-order phase transition that typically follows the nu- 
cleation and growth scenario [1]. While the stable crys- 
talline phase is favored energetically, the interface with 
the disordered phase penalizes small nuclei and leads to 
a free energy barrier. Hence, for crystallization to com- 
mence, one has to wait for a rare fluctuation that results 
in a large enough critical nucleus, which subsequently 
grows until it spans the system. This implies that by 
avoiding nucleation the density of a suspension can be 
increased beyond the thermodynamic freezing density. 
Such a metastable state is called supersaturated. Simple 
expressions for the free energy barrier can be obtained 
within classical nucleation theory (CNT) [2-4], which as- 
sumes a spherical nucleus. 

Subjecting a colloidal suspension to shear flow drives 
the system out of thermal equilibrium. Consequently, 
the concept of a free-energy is no longer well-defined 
and, strictly speaking, CNT cannot be applied anymore. 
Nevertheless, sheared colloidal suspensions still crystal- 
lize but with nucleation and growth kinetics that may 
be significantly altered from those in equilibrium [5, 6]. 
One might even find a dynamical coexistence of liquid 
and solid phases [7, 8]. Previous studies of the effects 
of shear flow on the crystallization rate are not con- 
clusive: On one side, shear-enhanced crystallization has 
been reported for experiments [9-13] and numerical sim- 
ulations [14, 15]. On the other side, a suppression of 
nucleation has been observed experimentally [16] and nu- 
merically [17, 18]. Others report an optimal strain rate 
for crystallization in supersaturated hard-sphere-like sus- 
pensions [19, 20] and protein solutions [21]. An optimal 
strain rate has also been reported in two dimensions for 
numerical simulations of Yukawa- type [22], Ising [23], and 
depletion-driven attractive [22] systems, and in a three- 
dimensional model glass [24]. 

Colloidal systems have the advantage that both the 



spatial and temporal evolution can be monitored directly 
in experiments via light scattering [25-27], or in real 
space via confocal microscopy [20, 28], see also Ref. 29 
for a review. Moreover, the interaction potential between 
colloidal particles can be tuned [30] . In the limit of per- 
fect hard-spheres, the crystallization process is driven 
solely by entropy and the phase diagram depends on 
the density only. For repulsive interactions, temperature 
plays a role; the density is, however, still the dominant 
control parameter. For a numerical treatment, molecu- 
lar dynamics as well as Brownian dynamics simulations 
appear to be ideally suited to track the temporal de- 
velopment of crystalline nuclei since particle coordinates 
are accessible and arbitrary interaction potentials can be 
used. However, since nucleation is hampered by a large 
free-energy barrier, it is a rare event. Depending on the 
height of the barrier, sampling such events can range from 
difficult to prohibitive unless specifically tailored meth- 
ods are employed. Such methods include transition path 
sampling [31, 32], forward flux sampling [33, 34], and 
umbrella sampling [35, 36]. 

In this work, we present numerical results for a model 
colloidal suspension in three dimensions, where parti- 
cles interact via hard-core exclusion plus a short-ranged 
Yukawa repulsion. The suspension is strongly super- 
saturated such that the probability for nucleation is 
shifted into a regime accessible by straightforward meth- 
ods without making use of importance sampling schemes. 
Here, we employ undcrdamped Langevin dynamics. Al- 
though the suspension is strongly supersaturated, we ob- 
serve that crystallization still proceeds via nucleation and 
growth. For the densities studied, small clusters appear 
and disappear until a cluster reaches the critical size. The 
following growth process is dominated by a single large 
cluster. We find that the time for this largest cluster to 
reach a specified size (larger than the critical size but 
smaller than the total number of particles) decreases for 
small strain rates but strongly increases at larger strain 
rates. There is, therefore, an optimal strain rate for the 
crystallization process. The reason is that shear flow is 
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responsible for two competing effects: suppression of the 
initial formation of a critical nucleus and enhancement 
of its growth once it has formed. Strain rates consid- 
ered here are so low that shear-induced layering plays no 
role [9, 37]. 

For unsheared liquids it has recently been emphasized 
that nucleation resembles a two-stage process and that 
pre-structuring of the liquid plays a crucial role in the 
formation of the critical nucleus [38-41] : Droplets of the 
stable phase appear preferentially in regions of the su- 
persaturated liquid that are still amorphous but where 
particle have already developed a loose connectivity with 
their neighbors. We confirm the two-stage scenario for 
the model studied here and demonstrate that the main ef- 
fect of the shear flow is to disrupt the formation of such a 
pre-structured liquid. Hence, we find that the suppressed 
nucleation under shear flow has its origin in the inhibition 
of structuring in the liquid rather than in the "demoli- 
tion" of crystalline clusters. The enhanced growth rate of 
sufficiently large clusters can be attributed to either con- 
vection or a faster reorganization of the cluster, making 
it thus easier to incorporate new particles. 

The paper is organized as follows. In Sec. II, we present 
the system and simulation details. Moreover, we describe 
the structural order parameters that we employ to deter- 
mine the structure of the local environment of single par- 
ticles. In Sec. Ill, we present and discuss our simulation 
results for both the nucleation and growth stage before 
we conclude in Sec. IV. 



II. METHODS 
A. System and simulation details 

We consider a mono-disperse colloidal suspension con- 
sisting of N = 4860 particles in a simulation box of con- 
stant volume V with volume fraction cj) = nN/ (6V) . We 
study two densities with corresponding box dimensions 
17.709x 17.351 x 18.404 for <j) = 0.45 and 16.876x 16.535 x 
17.538 for <f> = 52. Throughout the manuscript, we mea- 
sure length in units of the particle diameter a and energy 
in units of the thermal energy fceT, where T is the tem- 
perature of the suspension and &b is Boltzmann's con- 
stant. The time it takes for a particle to diffuse a distance 
corresponding to its diameter a defines the unit of time 
37t?7o 3 /(A;bT), where r\ is the viscosity of the solvent. 

We employ underdamped Langevin dynamics given by 
r fe = v fe and 



Thermal fluctuations are modeled by the stochastic forces 
£ with zero mean and correlations 



mvfc = -V k U - [v fe - u(r fe )] + £ k , 



(1) 



with positions and velocities Vfe. We choose the dimen- 
sionless mass m = 1 to be unity such that the relaxation 
time of the momenta equals the diffusive time-scale. The 
interaction forces are described by the potential 



U = J2u(\r i -r j \). 



(&(t)€j(0)) = 215 tf J(i), 



(3) 



where Sij denotes the Kronecker symbol, S(t) the Dirac 
distribution, and 1 the identity matrix. Moreover, we 
impose an external linear shear flow by means of a sol- 
vent velocity field u(r) = jye x entering the friction 
term, where 7 is the strain rate and e x the unit vec- 
tor in x direction. In the simulations, we employ peri- 
odic boundary conditions using the Lees-Edwards sliding 
bricks method [42]. 

Particles interact pairwise via the repulsive Yukawa 
potential 
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with hard-core exclusion. The strength of the potential 
is given by the energy at contact e, and its range is deter- 
mined by the inverse screening length k. The magnitude 
of k is mainly influenced by the ion concentration in the 
solvent and interpolates between Coulombic (low k, low 
ion concentration) and hard-sphere interactions (large n, 
high ion concentration). In this study, we choose e = 10 
and k = 8. For this set of parameters, the freezing vol- 
ume fraction is 4>* ~ 0.38. We are interested in the in- 
fluence of a weak but steady shear flow on the crystal- 
lization dynamics of colloidal suspensions under highly 
supersaturated conditions (</> = 0.45 and 0.52). 

The equations of motion (1) are integrated using a 
version of the velocity Verlet algorithm with time step 
At = 5 x 10~ 4 . Initial configurations for different runs 
are generated by equilibrating the system at low densi- 
ties (4> = 0.2). The volume fraction is then increased 
stepwise by rescaling the simulation box and the parti- 
cle coordinates until the final volume fraction is reached. 
We use 20000 time steps. Once we reach the final volume 
fraction, we switch on the shear flow with strain rate 7. 



B. Structure analysis 

In order to describe the crystallization process quan- 
titatively, we need a way to distinguish between liquid 
and solid structures. For the hard-core Yukawa system, 
the phase diagram of the equilibrium bulk structure in- 
cludes, beside the liquid phase, the two crystalline struc- 
tures body-centered cubic (bec) and face-centered cubic 
(fee) [1, 43, 44]. In the limit of hard-sphere interac- 
tion (high k), the free-energy difference between an fee 
and a hexagonal close-packed (hep) configuration is very 
low [45]. Therefore, hep structures are likely to occur 
along with fee structures, as has been observed in mi- 
crogravity experiments [46]. For the parameters stud- 
ied here, the bulk equilibrium structure is fee. However, 
in the spirit of the Ostwald step rule [47], intermedi- 
ate structures may be of a different type. Small nuclei, 
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e.g., have been found to belong predominantly to the bcc 
structure in a Lennard- Jones liquid [48]. Hence, in our 
analysis, we will not only discern fee from the liquid state, 
but also include hep and bec-structures. 

We employ different variants of the well-known Stein- 
hardt order parameters [49] to determine the local envi- 
ronment of a single particle. The basic idea is to con- 
struct quantities sensitive to the rotational symmetry of 
the local environment of the particles. To that end, one 
locates the set of neighbors Nb(k) of the fcth particle 
with size |iVb(fc)| and computes the connecting vectors 
Ykj = rfc — Yj of the central particle with its jth neighbor, 
where a neighbor is defined as another particle within a 
range not exceeding i?b- This range is frequently chosen 
as the minimum between the first and second shell in the 
pair correlation function. We define the complex vector 
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where the functions Yi m are spherical harmonics and I > 
0, |m| ^5 I. The vectors qi m depend sensitively on the 
choice of I. 

The first step is two distinguish fluid particles (disor- 
dered environment) from solid particles (ordered environ- 
ment). To this end, we make use of a recently introduced 
variant of the Steinhardt order parameters that averages 
over the second-neighbor shell [50], 
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where 
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where N^(k) is the set of neighboring particles including 
the fcth particle itself. Averaging the order parameter 
this way sharpens the distinction between different struc- 
tures at the expense of spatial resolution. For I — 6, the 
probability distributions of qg for fluid and solid particles 
are well separated, providing a good way to discriminate 
these two basic structure types from each other (data not 
shown). We regard a particle as fluid if q@ < 0.4 and as 
solid otherwise. 

Having determined these two particle sets, we further 
split the fluid particles into two subsets: liquid and pre- 
structured. While liquid particles have a truly disordered 
environment, we identify particles as pre-structured that 
have an environment that does not qualify as solid but 
where nevertheless some "bonds" between particles have 
formed. To concretize the concept of a bond, we consider 
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FIG. 1: Scatter plot of W4-we trajectories for perfect hep, bcc, 
and fee crystals subject to thermal fluctuations. The W4-wq 
plane is divided into the indicated regions, which are used to 
determine the crystal structure of solid particles. 
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with q* the complex conjugate of q. This product defines 
a measure for the strength of the correlation between 
the surrounding structures of the fcth and the jth parti- 
cle. We regard two neighboring particles as "bonded" if 
Sfcj > 0.5 [51] and denote the number of bonds for the 
fcth particle as nt>onds- While also in the liquid particles 
will have bonds, we consider particles that have at least 
8 bonds with neighboring particles (but still q & < 0.4) as 
pre-structured. 

Finally, by employing another type of averaged order 
parameters [50] 
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we are able to discern the different crystalline structures 
within the solid particles. The term in brackets is the 
Wigner-3-j symbol, which is related to Clebsch-Gordan 
coefficients. The sum runs over all combinations of mi, 
m2 and ?Ti3 with mi + m 2 + m 3 = 0. Using two param- 
eters Wi and u)6, we obtain a good separation between 
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FIG. 2: Decision tree for the assignment of a structure to 
a particle. On the first level, the particle is classified either 
as fluid or solid. On the second level, the fluid particles are 
split into pre-structured and liquid particles, while the solid 
particles are subdivided into the crystalline structures hep, 
bec, and fee. The double-headed arrows indicate the transi- 
tions we observe primarily during the crystallization process: 
While among the solid structures all transitions occur, the 
liquid particles almost never reach a crystalline state without 
passing through a pre-structured configuration first. 



the distributions in the w^wq plane, see scatter plot in 
Fig. 1. The id 4 distribution is widely separated for hep 
and fee structures, while w$ separates bec from hep and 
fee. Hence, a solid particle is classified as fee if W4 < 
and wq < —0.005 (lower left region), as hep if W4 > and 
wq < (lower right region), and as bec otherwise (upper 
region). Fig. 2 illustrates the complete classification pro- 
cess. Moreover, it shows the possible transitions between 
the different structures occurring during the crystalliza- 
tion. Between the solid states, each of the transitions 
occurs. However, the transitions are biased towards an 
fee environment. As we will see in the following, direct 
transitions from the liquid into a crystalline state are 
rare and the solidification of a liquid particle generally 
advances through a pre-structured state first. 



C. Discrete state model 

The structure analysis enables us to categorize sin- 
gle particles according to their local environment. For a 
given configuration of particle positions at time t, we de- 
fine an indicator function hk (t) for every particle k which 
takes on one of the five values: liquid (liq), pre-structured 
(pre), hep, bec, or fee. The population (fraction of par- 
ticles) of structure type i is 
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As the suspension evolves, the structural environment of 
particles will of course change. To quantify these changes 



we define the fluxes 
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which count the number of particles that have been con- 
verted from structure i into structure j within the time 
interval [t,t + St]. If not indicated otherwise, we set 
St = lOOAi. 

In the following, we consider the size n of the largest 
cluster as an order parameter characterizing the progress 
of the crystallization process. Clusters are identified as 
the sets of all solid particles that are mutually bonded 
(in the sense defined above that S$ > 0.5). We de- 
fine the average population at fixed cluster size n in the 
suspension as the conditional average 



d(n) = (ci(t)S nin ^) /Z n 



(12) 



where (•) averages over time and over different realiza- 
tions of the crystallization process. Hence, Z n = (<5„.„(t)) 
counts how many times the largest cluster size n occurs 
in all runs considered. Moreover, we define a 5 x 5 trans- 
fer matrix T, the components of which are given by the 
fluxes as 
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with normalization J^ . = 1. The component 

Ti^j{n) of this stochastic matrix quantifies the fraction 
of particles in state i that convert on average into state 
j within the time interval St. The eigenvalues of T can 
be sorted, Ao > Ai > ■ • • > A4, with Ao = 1. The compo- 
nents of the corresponding right-hand-side eigenvectors 
fulfill 
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The product T(n(t))c(t) = c(t + St) yields the average 
population a time St later with c(t) = (ci(t)). Apply- 
ing the transfer matrix repeatedly describes an effective 
Markovian dynamics at fixed cluster size. Under this dy- 
namics, the average population after a time r has elapsed 
becomes 



c(r) 



Ec« 



A 



(15) 



with implied time scales T a (n) = —St/]nX a (n) and coef- 
ficients ( a = w( Q ) • c(0). Hence, for r 3> n, the system 
approaches a (quasi)-stationary average population given 
by w^ ). The relaxation time is determined by t\. 

For a freely evolving suspension, the size of the largest 
cluster is of course not constrained to a fixed n. Sup- 
pose there is a time-scale separation: After a change 
of n the system relaxes into the quasi-stationary state 
before the cluster size changes again. Then the actual 
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FIG. 3: (a) Liquid, (b) pre-structured, and (c) solid fraction of 
the quasi-stationary populations wCI as function of n (lines) 
and the corresponding actual populations c (symbols). 



FIG. 4: (a) Average times of crystallization r x as a function 
of the shear rate with standard error bars for the estimation 
of the mean value, (b) Temporal development of n for typical 
runs at (f> = 0.52 for the strain rates 7 = 0.0, 0.02, and 0.08. 



average populations of local structures measured in the 
simulations should be equal to the stationary eigenvec- 
tor: c(n) ~ w'°'(n). Indeed, as shown in Fig. 3, this 
is the case to a very good degree. Shown are the ac- 
tual and quasi-stationary fractions of particles in a liq- 
uid, pre-structured, and solid local environment, where 

1 (0) (0) . (0) . (0) 

Csol = Chcp + Cbcc+Cf cc and w; o { = Whc p + W bcc + W fcc sum 

the contributions of all crystalline particles. Deviations 
are small but increase with increasing strain rate/cluster 
size up to maximal 7% at n = 1000. This demonstrates 
that the growth of the cluster is a slow process and that 
the lag time 5t is sufficient to sample the fast dynamics. 

III. RESULTS 

A. Crystallization rate 

We first examine the overall effect of shear flow on the 
total time it takes a supersaturated suspension to crystal- 
lize. According to the protocol described above, we start 
with a supersaturated, fluid suspension and apply shear 
flow at t = 0. Once the largest cluster in the suspension 
has reached a size of n — 800, we stop the simulation and 
denote the elapsed time as t x . The stopping size is clearly 
larger than the critical nucleus size but small enough 
(compared to the total number of particles) to minimize 
finite-size effects due to the periodic boundaries. Av- 



eraging over 500 independent runs for each strain rate, 
we obtain a non- monotonous dependence of t x on 7, see 
Fig. 4(a). For vanishing strain rate, the denser suspen- 
sion takes longer to crystallize. Increasing the strain rate, 
t x decreases for both densities. For the larger density, 
however, the decrease is much stronger and r x drops even 
below the value for the less dense system. In both cases, 
we observe an optimal strain rate 7* at which the crystal- 
lization process is fastest with 7* ~ 0.02 and 7* ~ 0.01 
for <f> — 0.52 and <f> = 0.45, respectively. At high strain 
rates, crystallization becomes rare. Already for 7 = 0.1 
and 7 = 0.05, for <f> = 0.52 and <f> = 0.45, respectively, a 
significant part of the 500 runs did not crystallize within 
2 x 10 6 time steps. We thus find an accelerated crystal- 
lization for small but nonvanishing strain rates, while for 
higher 7 crystallization is more and more suppressed. 

Following the development of the size n of the largest 
cluster in the suspension in time, we find qualitative dif- 
ferences in the way the crystalline state is reached for dif- 
ferent strain rates, see Fig. 4(b). Without shear flow, the 
system starts crystallizing quickly but the cluster grows 
slowly, whereas in the high-shear case, the system stays 
at a low n for some time before crystallization is initiated. 
Afterwards, the largest cluster develops rapidly. Close to 
the optimal strain rate, we find features of both limiting 
cases. The crystallization process starts almost as early 
as without shear flow, but progresses more rapidly later 
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FIG. 5: Average growth rate of the size n of the largest cluster 
at the nucleation stage (n small) as a function of n for different 
strain rates at (j) = 0.52. For large 7, the growth rate of the 
largest cluster is strongly suppressed in the interval 10 < n < 
40. 
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FIG. 6: Transfer matrix components for transitions from the 
liquid to the pre-structured state (solid lines) and the cor- 
responding backwards transitions (dashed lines) for volume 
fraction <j> = 0.52 as a function of cluster size n. 
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B. Shear flow suppresses nucleation 

1. Growth rate 

A pertinent quantity to study is the growth rate of the 
largest cluster 

»(n)= ^ t + 6t l- n ® 6 nMt) )/z n (16) 

averaged at fixed n, i. e., at a specific stage in the de- 
velopment of the largest cluster. The dot denotes the 
rate of change. In Fig. 5, we show this quantity for 
n < 50. For 7 = 0.02, the average growth rate is only 
marginally smaller than the one in the unshcarcd case. 
For 7 = 0.08, however, the growth is strongly suppressed 
in a broad range 10 < n < 40. This result could be due 
to a shear-induced inhibition of the formation, or due 
to a shear- induced destruction of small clusters. Either 
way, the likelihood for the formation of a critical cluster 
is strongly reduced. In the remainder of this subsection, 
we discuss the relevance of these mechanisms. 



2. Effect on pre-structured liquid 

Following the two-stage scenario for nucleation, a crys- 
talline cluster is not likely to occur in the middle of an 
entirely random distribution of particles. Rather, in a 
region of the liquid which has already acquired a loosely 
ordered state, fluctuations transforming parts of this pre- 
structured liquid into a crystal are much more likely to 
occur. This scenario also holds in the system studied 
here: While the fraction of pre-structured particles trans- 
ferring to a solid state is on the order of a few percent, 
7p ro -i.soi ~ 0.03, the corresponding fraction of particles 



converting directly from liquid to solid is smaller by more 
than three orders of magnitude, Tii q _j. so i < 10 -5 . 

In order to study the effect of shear flow on the struc- 
ture of the liquid, we record the transfer matrix compo- 
nents Xij q _,.p ro and corresponding backward component 
T pre _>iiq, see Fig. 6. Note that Ti-yj describes the fraction 
of particles in structure i transferring on average to struc- 
ture j in the following time interval St. Hence, the actual 
net current of particles changing from i to j depends as 
well on the population in these states. Hence, although 
31iq->pre < 2pre->Uqj the net current of particles is still 
directed from liquid to pre-structured, as the liquid state 
contains much more particles than the pre-structured one 
at this early stage in the crystallization process. In Fig. 6, 
we show Tjiq^pro and T pre ^ii q for three different strain 
rates. We find that shear flow has a significant influence 
on the development of structure in the liquid. On the one 
hand, the establishment of bonds is inhibited, as can be 
seen from the reduced values for Tii q _j. prG . On the other 
hand, structure in the liquid is destroyed, resulting in an 
enhanced value for T pre ^n q . Consequently, compared to 
the unsheared case, we find a much smaller fraction of 
pre-structured particles c pre both in the stationary and 
in the actual composition of the suspension for 7 = 0.08, 
see Fig. 3(b). In other words, shear flow prevents the 
liquid from developing a loose structure. 

In Fig. 7, we follow the evolution of the largest crys- 
talline cluster for strain rates 7 = 0.0 and 0.08. Clearly, 
in both cases the crystallization process is dominated 
by a single cluster. Note that crystalline clusters (large 
spheres in blue, gray, and red) are composed of different 
local structures, but there seems to be no tendency for a 
certain structure type to form a core or surface (e.g., a 
core of fee particles with a cluster-fluid interface formed 
by bec particles). Crystalline clusters are surrounded by 
pre-structured particle (small green spheres). Although 
most prominent in the vicinity of solid clusters, this loose 
structure can be found throughout the suspension. In the 
strongly sheared case, however, we observe considerably 
less pre-structured particles than in the unsheared sus- 
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FIG. 7: Snapshots of the suspension during the crystallization 
process for 7 = 0.0 (left column) and 7 = 0.08 (right column) 
for three different sizes of the largest cluster n. The large 
spheres are solid particles with hep (blue), bec (gray), and 
fec-structured (red) environments. Pre-structured particles 
are shown as small sphere (green), whereas liquid particles 
are not shown for clarity. 

pension, which indicates the shear-induced disruption of 
a loosely structured fluid. 

3. From pre-structured to solid 

Once a pre-structured but still amorphous environ- 
ment has formed, the pre-structured liquid has yet to 
transform into a crystalline cluster. Hence, the next step 
is to focus on the influence of shear flow on the second 
part of the transition from liquid to crystal. We trace 
the transitions between pre-structured and the crystalline 
structures hep, bec, and fee and show the corresponding 
transfer components in Fig. 8. Here, the influence 

of the shear flow is much smaller than for the transitions 
between liquid and pre-structured. The rates from the 
crystalline states to the pre-structured one describe the 
destruction of crystalline clusters. Interestingly, these 




100 100 100 200 



n n n 

FIG. 8: Transfer matrix components for transitions from the 
pre-structured state to hep, bec, and fee (from left to right, 
solid lines) scaled by a factor 10 and their corresponding back- 
wards fluxes (dashed) for 7 = 0.0, 0.02, and 0.08. 

are not enhanced by the shear flow but even somewhat 
reduced. Transitions into the different crystalline states 
are also barely affected by the shear flow. Hence, we 
find that the shear flow is not strong enough to actually 
destroy or shrink clusters once they have formed. 

C. Shear flow enhances growth of clusters 

1. Growth rate 

Complementary to Fig. 5, we now plot the growth 
rate Eq. (16) over a wider range of cluster sizes shown 
in Fig. 9. Note that for the larger clusters considered 
here, we need to increase the time interval over which 
the change in n is evaluated to St = 5000Ai in order to 
separate the growth trend from the fluctuations. Com- 
pared to the unsheared case, the growth rate is enhanced 
significantly once the shear flow is switched on. Fur- 
thermore, for cluster sizes n > 400 the growth rate is 
dominated by a linear term, the slopes of which them- 
selves grow proportional with the strain rate, see inset in 
Fig. 9. We observe the same behavior for <fr — 0.45 (data 
not shown) . The dominant contribution of the shear flow 
to the growth rate for n > 400 is thus a linear term of 
the form 

v[n) = Bjn, (17) 

where the proportionality constants for the slopes B can 
be determined from least-square fits to the data. We 
find B = 2.28 ± 0.10 and 2.75 ± 0.16 for <j) = 0.52 and 
4> = 0.45, respectively. This strongly shear-dependent 
growth rate reflects the enhanced cluster growth for large 
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FIG. 9: Average growth rate of the largest cluster as a func- 
tion of n for different strain rates for volume fraction <j> = 0.52. 
The solid lines are linear fits to the data for n > 400. Inset: 
Slope of the fit functions as a function of strain rate 7. The 
solid line is a linear fit starting at the origin. 



times as shown in Fig. 4(b) on a more systematic level. 
We now discuss two possible mechanisms of how shear 
flow influences the growth stage of crystallization. 



2. Mechanism I: Convection 



3. Mechanism II: bcc grows fastest 

A second explanation for an enhanced growth rate 
might be that one of the local structures can be grown 
faster. In Fig. 10(a)-(c), we show the components of the 
transfer matrix for transitions between the three solid 
structures. The most prominent effect of the shear flow is 
that the largest strain rate strongly facilitates transitions 
towards the bcc structure. For the other transitions, the 
effect of the shear flow is weaker. The rates towards the 
hep structure grow slightly with increasing 7 and the op- 
timal strain rate 7 = 0.02 enhances transitions from bcc 
to fee. 

Fig. 10(d) shows the relative composition of all solid 
particles in the suspension as a function of cluster size 
n. For small 7, most solid particles belong to fee, which 
is the stable equilibrium structure. Interestingly, even 
in the unsheared case does the fraction of fee particles 
slightly decrease as the cluster becomes bigger. In consis- 
tency with the shear dependence of the transition rates, 
at the optimal strain rate 7 = 0.02, the fraction of bcc 
particles is decreased and the number of fee particles is 
increased, whereas the shear flow has almost no influ- 
ence of the fraction of hep particles. Quite in contrast, 
the bcc particles overtake the fee particles at the higher 
strain rate 7 = 0.08 and become the dominant structural 
type. As can be seen in the snapshots in Fig. 7, shear 
flow allows for much larger bcc domains. This can be 
interpreted as a consequence of the Ostwald step rule: 
in the unsheared case bcc is more easily formed but fee 
is the more stable structure. Shear flow stabilizes bcc 
in relation to fee and hep, thus allowing for larger frac- 
tions of this structure type. If bcc is indeed the fastest 
growing structure, this mechanism would also lead to an 
enhanced crystalline growth. 



A growth rate that is a linear function of n might arise 
from convection, which can be understood as follows. We 
assume that particles in the vicinity of a crystalline nu- 
cleus are more likely to crystallize. Shear flow enhances 
the number of particles passing through this direct vicin- 
ity of the nucleus which we model as a sphere with radius 
R. The particle current entering this sphere caused by 
the shear flow reads 

If 4 3 

2 Js Pl ^ yGx ' dS ' = 3 PlR ^ ocA/n ( 18 ) 

with the surface of the sphere S and the number density 
of particles in the liquid pi = 7T0/6. The particle cur- 
rent is thus proportional to the size of the cluster, and 
proportional to the strain rate 7. Hence, we obtain the 
functional form of Eq. (17) where the free parameter B 
takes into account deviations from the spherical shapes 
of the cluster and the probability with which particles 
become solid (attach to the cluster). 



IV. CONCLUSIONS 

We have performed Langevin dynamics simulations to 
investigate the overall effect of weak linear shear flow on 
the crystallization process in a model colloidal suspen- 
sion. We found that the time it takes for the suspension 
to become solid exhibits a minimum at a finite strain 
rate. This can be explained as the result of two coun- 
teracting effects: at the early nucleation stage, the shear 
flow inhibits the formation of a critical nucleus, while at 
the later stage growth is enhanced. Both effects become 
manifest in the growth rate of the largest cluster, see 
Fig. 5 and Fig. 9. For larger strain rates the growth rate 
of small clusters is vanishingly small, whereas already 
small strain rates significantly speed up the growth pro- 
cess of the crystalline clusters after nucleation. 

To gain further insight into the different mechanism, 
we have employed a discrete state model. The state space 
of this discrete model comprises the relevant local struc- 
tural environments of a single particle. The microscopic 
environment of a particle might be fluid or solid. We 
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FIG. 10: Transfer matrix components for the transitions be- 
tween the three solid local structures: (a) hep and bec, (b) bec 
and fee, and (c) hep and fee. (d) Relative populations of solid 
particles in structure hep, bec, and fee Ci/c BO \ measured in the 
simulations. 



further distinguish between the unordered liquid and a 
pre-structured liquid, which still is amorphous but where 
already a loose structure between neighbors has formed 
as determined from a bond criterion. The solid particles 
are classified according to their crystal structure as either 
fee, bec, or hep; leading to a total of five discrete states. 
From the simulation data we can measure the popula- 
tions, i.e., the fraction of particles that resides in each 
structure. Since the simulations advance from the initial 
liquid towards the full crystal, these populations change 
over time. We have found it convenient to calculate con- 
ditional averages with respect to the size of the largest 
cluster. Moreover, we have determined the fluxes be- 
tween the discrete states, from which we built a transfer 
matrix. The quasi-stationary state of this transfer matrix 



agrees well with the measured actual populations. 

Nucleation in our model system clearly proceeds from 
liquid through pre-structured to one of the crystalline 
structures. The effect of shear flow is to disrupt the for- 
mation of the pre-structured regions, see Fig. 6. Conse- 
quently, there are fewer pre-structured particles and the 
probability to form a critical nucleus is drastically low- 
ered. On the other hand, once a critical droplet of solid 
particles has formed, the growth rate is enhanced by the 
shear flow. The emerging functional form (17) can be 
explained by convection alone: The flow field constantly 
changes the vicinity of crystalline clusters and is thus able 
to enhance the number of liquid particles under the "in- 
fluence" of a crystalline cluster. However, the shear flow 
also influences the composition of the cluster, hinting at 
a second mechanism. We found that shear flow facilitates 
transitions towards the bec structure and thus enhances 
the fraction of bec particles in the cluster. The stable 
bulk structure is fee and, therefore, there is a driving 
force towards fee. Destabilizing fee (or stabilizing bec) 
allows for larger domains of bec particles without the 
need to convert bec, which in turn allows the cluster to 
grow faster. 

Our data show that the depth of the minimum in the 
duration for the crystallization process is much smaller 
for a lower supersaturation. This point entails the ques- 
tion whether there is a lower bound on the density be- 
low which this minimum vanishes. How the existence 
and depth of this minimum depends on the supersatura- 
tion remains a topic for future investigation. Moreover, 
it would be interesting to study the influence of hydro- 
dynamic interactions on the crystallization process [52] 
under shear in general, and on the magnitude and the 
existence of such an optimal strain rate in particular. 
Beyond the insights we could gain for our specific model, 
we believe that the combination of techniques presented 
here will prove useful in the study of nucleation and re- 
lated processes. 
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